1 Workflow

  graph LR;
      M[OML]-->S[SPARQL];
      S-->API;
      API-->P[Python];
      P-->O[Oml2Owl];
      API-->MATLAB;
      MATLAB-->O;
      API-->MODELICA;
      MODELICA-->O;
      O-->M;

Below chart shows an actual workflow in this prototype.

flowchart LR
    id1(Build) --> id2(Query:SPARQL)
    id2 --> id3[massRollup:R]
    id3 --> id4[deltaV:Python]
    id4 --> id5[fuelMass:R]
    id5 --> id6(owl2oml)

2 First Prototype of Process and Data Flow

This prototype uses a Quarto notebook to orchestrate the different Like openMBEE, Jupyter or Quarto perform the execution of codes. Quarto has a function to execute Python and R together. We can bind the data from Python to R using the reticulate package.

flowchart LR
    id01[(OML)] --> id1
    id1(Build) --> id11[(owl)]
    id11 --> id2(Query:SPARQL)
    id2 --> id3[massRollup:R]
    id3 --> id4[deltaV:Python]
    id4 --> id5[fuelMass:R]
    id5 --> id6(owl2oml)
    id6 --> id01
    id2 -.-> id7((json))
    id7 -.-> id3
    id3 -.-> id8((df_mass))
    id8 -.-> id5
    id2 -.-> id9((orbit))
    id9 -.-> id4
    id4 -.-> id10((dv))
    id10 -.-> id5

In this approach, data is exchanged between different analysis tools like a rugby pass-and-receive.

2.1 load utility files

Code

library(reactable)
library(stringr)

searchDirectory <- function(iteration, pattern, parent_directory){
  for(i in 1:iteration){
      path <- list.files(parent_directory, recursive = TRUE, pattern = pattern, full.names = TRUE)
      if(length(path)){
        return(path)
      }
      parent_directory <- dirname(parent_directory)
  }
  print("file not found")
  return(path)
}

source(searchDirectory(4, "osr_common.R", (getwd())))
source(searchDirectory(4, "massRollupKepler16b.R", (getwd())))

Code

library(reticulate)
use_condaenv("py39")

2.2 Build and Run a Fuseki SPARQL endpoint for oml model query

Code

library(omlhashiR)
## oml_repository <- "../open-source-rover/"
oml_repository <- omlrepo
ret1 <- omlhashiR::oml_refresh()

Code

ret2 <- omlhashiR::oml_stop_Daemon(oml_repository)

Code

ret3 <- omlhashiR::oml_callTask(oml_repository, "build")

2.3 Start Fuseki Server

Code

ret4 <- omlhashiR::oml_callTask(oml_repository, "startFuseki")

2.4 Load in OWL ontologies to Fuseki Server

Code

ret5 <- omlhashiR::oml_callTask(oml_repository, "owlLoad")

2.5 Get Configuration & Mass Properties

Code

library(tansakusuR)
endpoint_url <- "http://localhost:3030/tutorial2/sparql"

configuration_root_iri <- "http://example.com/tutorial2/description/components#orbiter-spacecraft"

2.5.1 SPARQL CODE

Code

query_string <- '
PREFIX mission:     <http://imce.jpl.nasa.gov/foundation/mission#>
PREFIX base:        <http://imce.jpl.nasa.gov/foundation/base#>
PREFIX project:     <http://imce.jpl.nasa.gov/foundation/project#>
PREFIX vim4:        <http://bipm.org/jcgm/vim4#>

SELECT DISTINCT ?c1 ?c1_instancename ?c1_id ?c1_name ?c1_mass ?c2 ?c2_instancename ?c2_id ?c2_name ?c2_mass ?c3 ?c3_instancename ?c3_id ?c3_name ?c4
WHERE {

  VALUES ?c1 { <$configuration_root_iri> }


   OPTIONAL{ 
    ?c1 base:hasIdentifier ?c1_id ;
      base:hasCanonicalName ?c1_name ;
  }
    OPTIONAL {
        ?c1_mass_mag vim4:characterizes ?c1 ;
            vim4:hasDoubleNumber ?c1_mass .
    } 
  OPTIONAL{
    ?c1 base:contains ?c2 ;
    OPTIONAL{
      ?c2 base:hasIdentifier ?c2_id ;
          base:hasCanonicalName ?c2_name .
    }
    OPTIONAL {
        ?c2_mass_mag vim4:characterizes ?c2 ;
            vim4:hasDoubleNumber ?c2_mass .
    }     
    OPTIONAL{
      ?c2 project:isSuppliedBy ?c3 ;
          OPTIONAL{
            ?c3 base:hasIdentifier ?c3_id ;
                base:hasCanonicalName ?c3_name;
                project:isAuthorizedBy ?c4 .
          }          
    }    
  }


  BIND(STRAFTER(STR(?c1), "#") AS ?c1_instancename) .
  BIND(STRAFTER(STR(?c2), "#") AS ?c2_instancename) .
  BIND(STRAFTER(STR(?c3), "#") AS ?c3_instancename) .
 }
ORDER BY ?c2_id

'
query_string <- str_replace(query_string, "\\$configuration_root_iri", configuration_root_iri)

2.5.2 Query

Code

df_query <- send_query(endpoint_url,query_string)

2.6 Binding Variables: toR Parameter Table

OML descriptions shall specify this parameter table using IMCE Analysis Vocabularies.
OML descriptions shall specify mapping of binding variables for each simulation.
Get initial value from oml descriptions by SPARQL Query.

2.6.1 SPARQL CODE

Code

query_param_string <- '
PREFIX base:        <http://imce.jpl.nasa.gov/foundation/base#>
PREFIX mission:        <http://imce.jpl.nasa.gov/foundation/mission#>
PREFIX analysis:        <http://imce.jpl.nasa.gov/foundation/analysis#>
PREFIX sa:        <http://example.com/tutorial2/vocabulary/stateanalysis#>

SELECT DISTINCT ?iri ?statevariable ?value
WHERE {

#  VALUES ?analysisTarget {<http://example.com/tutorial2/description/statedictionary#orbiter-spacecraft.delta-v>}
   VALUES ?analysisTarget {<http://example.com/tutorial2/description/statedictionary#orbiter-spacecraft.mass.wet>}
  VALUES ?componentType { sa:StateVariable }
  VALUES ?relationType { sa:isAffectedBy }

  # Recursive Query to get state variables
  ?analysisTarget sa:isAffectedBy* ?iri ;
#  ?iri a ?componentType
  
  OPTIONAL{
    ?iri sa:hasStateValue ?value;
         }



  BIND(STRAFTER(STR(?iri), "#") AS ?statevariable) .
  BIND(STRAFTER(STR(?targetIri), "#") AS ?target) .
  BIND(STRAFTER(STR(?cpIri), "#") AS ?component) .
 }
ORDER BY ?iri

'

2.6.2 Query

Code

df_query_param <- send_query(endpoint_url,query_param_string)

2.6.3 Binding to R Dataframe

Code

# mapping to parameters to instance
# df_parameters <- data.frame(
#   initOrbit = as.double(df_query_param$value[df_query_param$statevariable=="orbiter-spacecraft.initial-orbit"]),
#   targetOrbit = as.double(df_query_param$value[df_query_param$statevariable=="orbiter-spacecraft.target-orbit"]),
#   m_wet = 0.0,
#   m_dry = 0.0,
#   m_fuel = 0.0,
#   dv = 0.0,
#   I_sp = as.double(df_query_param$value[df_query_param$statevariable=="orbiter-spacecraft.isp"])
# )

mapping <- list(c("dv", "orbiter-spacecraft.delta-v"),
  c("initOrbit", "orbiter-spacecraft.initial-orbit"),
  c("I_sp","orbiter-spacecraft.isp"),
  c("m_dry","orbiter-spacecraft.mass.dry"),     
  c("m_fuel","orbiter-spacecraft.mass.fuel"),
  c("m_wet","orbiter-spacecraft.mass.wet"),
  c("targetOrbit","orbiter-spacecraft.target-orbit")
)

# Convert to dataframe
mapping_df <- mapping %>% 
  map_df(~data.frame(.x[1], .x[2])) 
colnames(mapping_df) <- c("parameter", "statevariable")

df_parameters_before <- left_join(df_query_param, mapping_df, by=c("statevariable"="statevariable")) 

# parameter database for quarto with R
df_parameters <- df_parameters_before %>% 
  mutate(value2 = as.double(replace(value, is.na(value), 0))) %>%
  select(parameter,value2) %>%
  pivot_wider(
#    id_cols = age_cat,
    names_from = parameter,
    values_from = value2
    )

Code

df_table <- df_parameters %>% 
  mutate(id=1) %>%
  pivot_longer(!id, names_to = "parameter", values_to = "value") %>%
  select(-id)

datatable(df_table, options = list(autoWidth = TRUE, pageLength = -1))

2.7 MassRollUp

2.7.1 Tidy Data

Code

df_query$c1_type <- "component"

df_query <- df_query %>%
  mutate(c2_label = paste0(c2_id,": ",c2_instancename)) %>%
  mutate(c1_label = paste0(c1_id,": ",c1_instancename)) 

df_config <- df_query

# just for vis add NA to root
df_vis <- df_config %>%
  add_row(c1=NA,
          c1_instancename=NA,
          c1_id=NA,
          c1_name=NA,
          c1_mass=NA,
          c2=df_config$c1[1],
          c2_instancename=df_config$c1_instancename[1],
          c2_id=df_config$c1_id[1],
          c2_name=df_config$c1_name[1],
          c2_mass=df_config$c1_mass[1],
          c1_type=df_config$c1_type[1],
          c2_label=NA,
          c1_label=df_config$c1_label[1],  
          .before = 1)

plotCollapsibleTreeFromDataframe(df_vis, palette="BluYl", 
                                 parent="c2_label",
                                 child="c1_label",
                                 type="c2_id")

2.7.2 Prepare a Graph Data for mass rollup calculation

Code

library(igraph)

df_g <- df_config %>%
  mutate(parent = c1_instancename) %>%
  mutate(child = c2_instancename) %>%
  mutate(type = c2_instancename) %>%
  # mutate(parent = paste0(p_id, ": ", p_instancename)) %>%
  # mutate(child = paste0(c_id, ": ", c_instancename)) %>%  
  select("parent","child","type")

# df_g <- df_g[-1,]

g <- graph_from_data_frame(df_g, 
                           directed = TRUE, 
                           vertices = NULL)

2.7.3 RollUp Pattern

Code

root <- V(g)[1]

# Depth-first search is an algorithm to traverse a graph. It starts from a root vertex and tries to go quickly as far from as possible.
order <- dfs(g, V(g)[root], order.out = TRUE)$order

df_mel_before <- igraph::as_data_frame(g, what = "vertices") %>%
  arrange(factor(name, levels = names(order)))%>%
  filter(name !="NA")

# ここに、df_mass として massの値を入れたコンフィグを作る必要がある。

df_mass <-  df_config

df_mass <-  df_mass %>%
  mutate(mass = as.numeric(c2_mass)) %>%
  select(c2,
         c2_instancename,
         c2_id,
         c2_name,
         c1_id,
         mass) %>%  
  add_row(c2=df_mass$c1[1],
          c2_instancename=df_mass$c1_instancename[1],
          c2_id=df_mass$c1_id[1],
          c2_name=df_mass$c1_name[1],
          c1_id=NA,
          mass=as.numeric(df_mass$c1_mass[1]),
          .before = 1)%>%
  arrange(c2_id)


colnames(df_mass) <- c("c_iri","c_instancename","c_id","c_name", "p_id", "mass")

df_mel_before <- left_join(df_mel_before, df_mass, by = c("name"="c_instancename"))

Code

namekey="name"
masskey="mass"


df_mass_update <- massRollUp(g, root, df_mel_before, namekey = "name",masskey = "mass") %>%
  select(name, mass)


# リーフ数を取得
df_deg <- data.frame(
  name = names(degree(g)),
  degree = degree(g)-1, # num of components is degree(g)-1
  distance = distances(g)[V(g)[1],]
)

# distanceの値で、componentか、workpackageかを切り分けることができる

df_mel_after <- df_mel_before %>% select(-mass)
df_mel_after <- left_join(df_mel_after, df_mass_update, by=c("name"="name")) 
df_mel_after <- left_join(df_mel_after, df_deg, by=c("name"="name")) 

df_mel_after$componenttype <- ifelse(df_mel_after$distance == 0, "system", 
                            ifelse(df_mel_after$distance == 1, "subsystem", 
                                   ifelse(df_mel_after$distance == 2, "assembly", NA)))

Code

df_table <- df_mel_after %>%
  mutate(totalmass=mass) %>%
  select(c_id, name, totalmass, componenttype) 

df_table$componenttype <- factor(df_table$componenttype)

datatable(df_table, options = list(autoWidth = FALSE, pageLength = -1), filter = list(
  position = 'top', clear = FALSE
))

2.8 Binding Variables : toR

Code

df_parameters$m_dry <- df_mel_after$mass[1]
# df_parameters$m_dry <- 325.47880079

2.9 Binding Variables: toPython

Code

orbit_init=r.df_parameters.initOrbit.values.tolist()[0]
orbit_target=r.df_parameters.targetOrbit.values.tolist()[0]
m_final=r.df_parameters.m_dry.values.tolist()[0]
I_sp=r.df_parameters.I_sp.values.tolist()[0]
# orbit_init=r.initOrbit
# orbit_target=r.targetOrbit
# m_init=1000
# I_sp=r.I_sp

2.10 deltaV ( Using Python)

Code

import os
import sys

print('getcwd:      ', os.getcwd())

getcwd:       /home/runner/work/kepler16b-using-imce-vocabulary/kepler16b-using-imce-vocabulary/src/analysis

Code

sys.path.append(os.getcwd())
sys.path.append(os.getcwd()+'/src/analysis')

# from src.analysis import analysisOrbit
import analysisOrbit
from astropy import units as u

analysisRocket = analysisOrbit.hohmanTransfer(orbit_init=orbit_init, orbit_target=orbit_target)
# analysisRocket = analysisOrbit.hohmanTransfer(orbit_init=r.initOrbit, orbit_target=r.targetOrbit, m_init=r.m_dry, I_sp=r.I_sp)
# analysisRocket = analysisOrbit.hohmanTransfer(orbit_init=400, orbit_target=35786, m_init=5000, I_sp=300)

total_delta_v = analysisRocket.calculate_delta_v()


print(f"Total delta-v: {total_delta_v}")

Total delta-v: 3853.959363102248 m / s

2.11 Binding Variables: toR

Code

df_parameters$dv <- py$total_delta_v /1000.0 # km/s

2.12 fuelMass (Using R)

Code

source(searchDirectory(4, "calcWetMass.R", (getwd())))

# I_sp <- 350
# m_dry <- 325.47880079
# df_parameters$dv <- sum(df_parameters$dv)
# Inverse Problem
df_parameters$m_wet <- calcWetMass(df_parameters$dv, df_parameters$m_dry, df_parameters$I_sp)
df_parameters$m_fuel <- df_parameters$m_wet - df_parameters$m_dry
c(df_parameters$m_wet, df_parameters$m_dry, df_parameters$m_fuel)

[1] 6012.68 1957.00 4055.68

2.13 Finally, Update OML

2.13.1 parameter table updated

Code

df_table <- df_parameters %>% 
  mutate(id=1) %>%
  pivot_longer(!id, names_to = "parameter", values_to = "value") %>%
  select(-id)

datatable(df_table, options = list(autoWidth = TRUE, pageLength = -1))

2.13.2 parameters before and after

Code

df_parameters_after <-  left_join(df_parameters_before, df_table,
                                  by = c("parameter"="parameter"),
                                  suffix=c("_before","_after")) %>%
  mutate(status = case_when(
    as.double(value_before) == value_after ~ "unchanged",
    TRUE ~ "changed"
  ))

2.13.3 sending update query

Code

source(searchDirectory(4, "updateSparqlQuery.R", (getwd())))

df_update <- df_parameters_after %>% 
  filter(status == "changed")

endpoint_url <- "http://localhost:3030/tutorial2-tdb/"

for( i in 1:nrow(df_update) ){
  df <- df_update[i,]
  update_value <- as.character(df$value_after)
  update_iri <- df$iri

  df_ret <- updateSparqlQuery(endpoint_url, update_value, update_iri)
  
  ret <- send_update(endpoint_url = endpoint_url, df_ret$query_string_delete_tdb)
  ret <- send_update(endpoint_url = endpoint_url, df_ret$query_string_insert_tdb)
  
}

# cat(df_ret$query_string_delete_tdb)
# cat(df_ret$query_string_insert_tdb)

2.14 Run Save in OML-Vision

Using OML-Vision GUI, run the save task. Then you can see the oml file is updated.

2.15 Verify Update

This can be verified when you re-load the load task and send a query of query_param_string.

Code

endpoint_url <- "http://localhost:3030/tutorial2/"
df_query <- send_query(endpoint_url, query_param_string)


datatable(df_query, options = list(autoWidth = FALSE, pageLength = -1))

3 Discussion: Proposed Process and Data Flow

Instead of using the passing and recieving approach, alternative way is using the owl fuseki database as a hub to exchange the data. OML descriptions specify the process and mapping of data exchanged between analyssis containers.

To do this, we need to specify the analysis container with

input
output
mapping of input and output with other containers

flowchart LR
    oml[(OML)] -.-> id1
    id1(Build) -.-> owl[(owl)]
    id1 --> id2
    owl -.-> id2(Query:SPARQL)
    id2 --> id3[massRollup:R]
    id3 -.-> owl
    id3 --> id4[deltaV:Python]
    id4 -.-> owl
    id4 --> id5[fuelMass:R]
    id5 -.-> owl
    id5 --> id6(owl2oml)
    owl -.-> id6
    id6 -.-> oml
    classDef codpy fill:#ff0000;
    classDef coder fill:#00ffff;
    class id3,id5 coder;
    class id4 codepy;

4 TODO

df_parameters is not good idea to bind the parameters I think this could be better.

Code

df_parameters2 <- data.frame(
  parameter = c("initOrbit","targetOrbit","I_sp","m_dry","m_fuel","m_wet","dv"),
  value = c(400,35786,350,
          df_parameters$m_dry, 
          df_parameters$m_fuel, 
          df_parameters$m_wet,
          df_parameters$dv),
  type = c("input","input","input","output","input","output","output")
)

# how to access data
df_parameters2$value[df_parameters2$parameter == "m_dry"]

[1] 1957

Here is how to access values. when we add the mapping descriptions to analysis codes, we might automate the process of data exchange.

Code

df = r.df_parameters2
df[df['parameter']=='I_sp'].value.values[0]

350.0

4.1 OPTION: m_fuel ( Using Python codes)

Code

analysisRocket.calculate_init_mass(I_sp, m_final)
print(f"Total delta-v: {total_delta_v}")

Total delta-v: 3853.959363102248 m / s

Code

print(f"Required fuel mass: {analysisRocket.m_fuel}\n final: {analysisRocket.m_final}\n init: {analysisRocket.m_init}")

Required fuel mass: 4055.680381236295
 final: 1957.0
 init: 6012.680381236295

4.2 OPTION: deltaV ( Using R codes)

Code

source(searchDirectory(4, "calcDeltaV.R", (getwd())))

dv <- calcDeltaV(initOrbit=df_parameters$initOrbit, targetOrbit = df_parameters$targetOrbit)

[1] 2.397473 1.456487